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Abstract 

We present an analysis of the existing data on charmonium hadro-production based on non- 
relativistic QCD (NRQCD) calculations at the next-to-leading order (NLO). All the data on J/tp and 
■0(25) production in fixed-target experiments and on pp collisions at low energy are included. We find 
that the amount of color octet contribution needed to describe the data is about 1/10 of that found 
at the Tevatron. 

1 Introduction 

The production of charmonium and bottomonium states in high-energy collisions has always been the 
subject of considerable interest 1 . From the experimental point of view, charmonium decays into lepton 
pairs offer very clean signatures that are used not only for triggering and calibration but also to perform 
important physics studies. The decay B — > J /ip + X, for instance, provides an easy handle to QCD 
studies of b production and to the precise determination of some of the CKM parameters, such as sin 2/3. 

The importance of quarkonium is widely recognized also by the theoretical community. Charmonium 
and bottomonium states offer a unique laboratory for testing our understanding of QCD, and in particular 
of the interplay between the perturbative and non-perturbative regimes, which describe the physics 
of heavy-quark creation and that of bound state formation, respectively. The so-called color-singlet 
model [3] has been superseded by a rigorous framework, based on the use of non-relativistic QCD 
(NRQCD) 4 , an effective field theory that consistently includes relativistic corrections and provides a 
solid ground for accurate theoretical analyses. 

1 Supported by the Bundesministerium fur Bildung und Forschung, FRG, under contract numbers 05-7MP25I. 
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However, despite the theoretical developments and successes, not all the predictions of the NRQCD 
factorization approach have been firmly established. The first example is the universality of the non- 
pertubative matrix elements, on which the predictive power of this approach relies. Measurements at the 
Tevatron in proton-proton collisions suggest larger values for the color-octet matrix elements than those 
obtained at HERA in electron-proton collisions 0. Even more problematic is the measurement of J ftp 
polarization at the Tevatron. NRQCD predicts a sizable transverse polarization for J ftp's at high-pr, 
in contrast with the latest data that now clearly indicate that J ftp's are not transversely polarized [5]. 
While there is no quantitative and accepted explanation for this behaviour, it is generally argued that 
because the quarkonium mass is still not very large with respect to the QCD scale, in particular for the 
charmonium system, non-factorizable corrections may not be suppressed enough and/or the expansions 
in NRQCD may not converge very well. 

In the light of the present uncertain status, detailed studies on the range of applicability of the 
NRQCD approach, above all for charmonium states, are certainly welcome. In this work, we perform 
the NRQCD analysis of charmonium production data as obtained from experiments at fixed-target 
(with the exclusion of that induced by pion beams) and from pp collisions at low energy. Our purpose is 
twofold. First, we present an up-to-date collection of the experimental data on charmonium production 
in fixed-target experiments. Second, we study whether data are consistent with the NRQCD approach 
and in particular we extract information on the color-octet contributions, to be compared with that 
obtained from other experiments. Our analysis represents an improvement with respect to the previous 
studies presented in Refs. [7J|S], both in terms of accuracy of the theoretical predictions and for the more 
complete and comprehensive treatment of the experimental data. The outline of the paper is as follows. 
In the next section, we briefly review the framework of the NRQCD approach and state the theoretical 
results and assumptions that enter into our predictions. In Section |21 we discuss the experimental data. 
In Section we describe the fit strategy and present the results. We draw our conclusions in the last 
section. 

2 The NRQCD approach 

In the NRQCD approach, the cross section for the production of a quarkonium state H in a nucleon- 
nucleon interaction is expressed as a sum of terms, each of which factors into a short-distance coefficient 
and a long-distance matrix element: 

a(pp^H + X) = Y J I dx 1 dx 2 f i/p f j/p £ a(ij -» QQ [n] + x) (O h [n]>, (1) 

where the indexes i,j run over all the partonic species and n denotes the color, spin and angular 
momentum state of an intermediate QQ pair. The short-distance cross section a can be calculated as a 
perturbative expansion in the strong coupling a s . The NRQCD matrix elements (0 H [n]} (see Ref. [3] 
for their definition) are related to the non-perturbative transition probabilities from the QQ state n 
into the quarkonium H. They scale according to a definite power of the intrinsic heavy-quark velocity v 
(v 2 ~ 0.3 for charmonium and v 2 ~ 0.1 for bottomonium) [jjj. The general expression Q is thus a double 
expansion in powers of a s and v. While a formal and general proof of Eq. (pQ) is still lacking, it has been 
recently shown ^U] that it holds for high-p^ quarkonium production up to its two-loop description. In 
this work, we simply assume that soft effects do not spoil factorization and Eq. Q holds true also for 
total cross sections. 

The color-singlet short distance coefficients for spin singlet S-wave ( 1 5'q 1 ^), P-waves and all the leading 
color-octet coefficients are known at NLO for both photon-proton and proton-proton collisions [111 
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H 


(Of) 






• = (O 8 [ 3 P {8) ])/m 2 


J/1> 


1.16 GeV 3 


1.19 • 10~ 2 GeV 3 


1.0- 


10~ 2 GeV 3 




0.76 GeV 3 


0.50 • 10- 2 GeV 3 


0.42 


• 10~ 2 GeV 3 


XdO 


0.11 GeV 5 


0.31 • 10~ 2 GeV 3 







Table 1: Reference NRQCD matrix elements for charmonium production. The color-singlet matrix 
elements are taken from the potential model calculation of |141 115j . The color-octet matrix elements 
have been extracted from the CDF data in Ref. |17j . 



I12| . The color-singlet coefficient for 3 S[ is known at NLO only for photon-proton collisions [El • ^ n 
this respect our analysis cannot be considered as fully at NLO and should be updated once the NLO 
calculation for the color singlet term will be available. On the other hand, this is not an important 
limitation to our results, as will be made clear in the following. 

The non-perturbative matrix elements have to be extracted from the data. For the color-singlet 
terms this is straightforward. It can be easily shown that, up to relativistic corrections of order t> 4 , they 
can be related to those appearing in the corresponding decay rates and therefore can be extracted from 
measurements of decay widths. On the other hand, color-octet matrix elements can only be extracted 
from production processes, such as photoproduction or hadroproduction. The factorization hypothesis 
implies that the values extracted from different experiments should be universal. Table ^ shows the 
results of a fit performed on the CDF charmonium data 16 , providing the leading color octet terms 
for J/tp, ip(2S) and Xcj production at the Tevatron energy [E|. In the case of S- waves, the fact that 
transverse momentum distributions coming from CP-even states ( 1 Sq 8 ^ and 3 Pj 8 ') and the 3 S^ have 
different shapes, has been exploited to obtain information on their relative size. An equally detailed 
information cannot be extracted from fixed-target results, normally limited to the total production 
rates. 

The analysis performed here is based on a code implementing the NLO calculations of Ref. |12j . 
For the theoretical inputs, we make the following choices. We use /Xq = 2m c with m c = 1.5 GeV as our 
central value for the renormalization, factorization and NRQCD scales (all taken equal). We estimate the 
associated uncertainty by varying the scales between /j,q and 4^o- The strong coupling constant as(mz) 
is tuned to the one used in the PDF sets, i.e., CTEQ6m [IB] and MRST2002nlo 19 . We exploit spin 
symmetry to reduce the number of independent non-perturbative matrix elements, 

(oiepj)) = (2j+i)(o 8 ^( 3 p )), 
(o^esi)) = (2j+i)(or , ( 3 s 1 )) , (2) 

(0^( 3 Pj)) = (2J + 1) (Of°( 3 P )) 

and consider only leading color octet corrections. We take the non-perturbative matrix elements collected 
in Tab. H as our reference values. 2 Color-singlet matrix elements are kept fixed. For the color-octet 
matrix elements we adopt the relative normalization as that obtained from Tevatron measurements, 

2 See also Ref. [20 for a more recent analysis. 
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but for the S-wave color-octet matrix-element two overall multiplicative numbers Aj/^ and Xm2S) are 
introduced to be fitted with the fixed-target data. The term (0g cO ( 3 Si)) is left fixed to its reference 
value extracted from the Tevatron data. With the above assumptions, the data are fitted with only 
two free parameters, the A's, which can be interpreted as the fractions of the "overall Tevatron octet 
contribution" for J ftp and ip(2S) necessary to explain the fixed-target data. 
The following theoretical expressions for the cross sections are used: 

0>(2S) = cr° (2S) , 

2 

<rj/f = + E Br (XcJ - J/M) <% cJ + Br (V(2S) -> J/^X) a^ (2S) , (3) 
j=o 

^(2g) 

tup — > 

where the superscript D refers to the direct contribution. For the branching ratios we use [2^ 



Br(?/>(2S) - 


» J/1>X) 


= 57.6% 


Br(x c o - 


- J/M) 


= 1.18% 


Br(x c i - 




= 31.6% 


Br(x C 2 - 


- J/M) 


= 20.2% 



(4) 



3 Present experimental situation 

The measurements of J/ip hadroproduction have been performed in a time period spanning about thirty 
years. Over such a long period, several different experimental techniques have been used and different 
input information was available at the time of the measurements. Therefore, comparing results of 
different experiments on an equal footing requires an update of the published numbers on several aspects. 
For example, the charmonium branching ratios have changed with time and the treatment of the nuclear 
effects are not homogeneous. In our compilations, we update all the measurements with the current best 
knowledge of branching ratios and nuclear effects. 

3.1 Compilation of J/ip cross sections 

The cross section for J/ip production on a nuclear target of mass number A is parametrized as 

° P j% = *J/0-A a . (5) 

For this comparison the most precise measurement of a = 0.96 ± 0.01 22 at xf(J/iP) — is used, 
assuming its independence of the cms-energy. If an experiment has published cross sections obtained 
from different targets, Eq. © is applied to obtain a combined result. This is the case for the experiments 
CERN-PS HI] (H, C, W), NA50 [SZI (Be, Al, Cu, Ag and W) and HERA-5 [Ml (C, Ti and W). As 
suggested in Ref. 36], we combine the results of NA51 [SSI (H, D targets) and NA38 [SO] (C, Al, Cu, 
W). 

The results collected in Tab.Elare updated with the latest branching fractions |21| (5.88 ±0.10% for 
J/ip — ► fj- + fi~ and 5.93 ± 0.10% for J/ip — ► e + e~). If only the forward cross section (xp(J/tp) > 0) is 
given, we assume a symmetric daj/^/dxF distribution and multiply the forward cross section by two to 
obtain the total cross section. The cross sections of NA38, NA50 and NA51 are only quoted for a limited 
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phase-space. In order to obtain total cross sections, the extrapolation discussed in the NA51 paper is 
performed [22|. 

Unless it is clearly stated in the publications, we assume that the quoted uncertainties on branching 
ratio and atomic mass dependence were taken into account in the systematic uncertainty. Therefore we 
consistently update them in Tab. [21 

The updated results for mid-rapidity and total cross sections are displayed in Fig. 
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Figure 1: Production cross sections in proton-induced interactions from Tab. 21 as a function of the cms- 
energy. Proton-proton (pp) and proton-nucleus (pA) measurements are indicated by different symbols. 
Left: differential cross sections da p N /dy at y=0; right: total cross sections for J/tp and tp(2S) production. 



3.2 Compilation of ip(2S) cross sections 

The procedure described in the previous section is applied to the published ip(2S) cross sections us- 
ing the most recent branching fractions [2^ (0.73 db 0.08% for tp(2S) — ► [i + [i~ and 0.755 db 0.031% for 
i/j(2S) e+e") and a = 0.934 ± 0.010 [22] at x F (ip{2S)) ~ 0. 

The results on the ip(2S) cross sections (see Fig. ^) and the ratios between ip(2S) and J/ip cross 
sections are listed in Tab. El 

3.3 Comments on the available data 

The J/V> cross sections have been usually measured on large samples by triggering on dilepton decays. 
Therefore the most precise data have uncertainties dominated by systematic errors. The measurements 
show a good overall consistency, even though some of the results are hardly compatible. For instance, 
the two results at 20.6 GeV (E331 |2H1 and E444 [201) differ h Y roughly 2a. The E705 M result at 23.8 
GeV exceeds the UA6 32 one at 24.3 GeV by more than 3<r, in contrast with the expectation (and the 
general trend) that the cross section should increase with the cms-energy. This is an indication that 
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Experiment 


Reaction 


Vs 


dy \y=° 












(GeV) 


(nb / nucleon) 


(nb/nucleon) 


CERN-PS 23 


P A 


6.8 




0.65 ± 0.06 


WA39 24 


PP 


8.7 




2.4 ± 1.2 


IHEP 25 


pBe 


11.5 


16 ± 5.2 


20 ± 5.2 


E331 26 


pBe 


16.8 


84 ± 20 


122 ± 40 


NA3 


27 




pPt 


16.8 




80 ± 15 


NA3 27 


pPt 


19.4 




110 ± 21 


NA3 


27 




PP 


19.4 




124 ± 22 


E331 


2* 




P C 


20.6 




256 ± 30 


E444 


2^ 




P C 


20.6 




166 ± 23 


ISR | 


30. 




PP 


23.0 


100 ± 77 




E705 


31 




pLi 


23.8 




267 ± 30 


UA6 


32 




PP 


24.3 


104 ± 19 


152 ± 20 


E288 




, 


pBe 


27.4 


131 ± 33 


204 ± 51 


E595 






pFe 


27.4 


187 ± 12 


306 ± 18 


NA38/51 35, 36 


pA 


29.1 


169 ± 13 


292 ± 64 


NA50 37 


pA 


29.1 


188 ± 14 


325 ± 67 


ISR | 


38 




PP 


30 


154 ± 42 




ISR | 


39 




PP 


30.6 


111 ± 30 




ISR 


30 




PP 


31 


142 ± 93 




E672/706 40 


pBe 


31.6 




274 ± 60 


E771 


41 




pSi 


38.8 


202 ± 17 


333 ± 25 


E789 42 


pAu 


38.8 


170 ± 30 


327 ± 56 


HERA-B 43 


pA 


41.6 


392 ± 51 


663 ± 87 


ISR 44 


PP 


52 


204 ± 85 


716 ± 303 


ISR | 






PP 


52 


216 ± 54 




ISR | 


39 




PP 


52.4 


185 ± 12 




ISR | 


38. 




PP 


53 


229 ± 52 




ISR 


30 




PP 


53 


280 ± 161 




ISR 


39 




PP 


62.7 


172 ± 15 




ISR 


30 




PP 


63 


538 ± 346 




ISR 


38 




PP 


63 


250 ± 56 




PHENIX 46 


PP 


200 


1051± 255 


4000 ± 938 



Table 2: Updated differential {daji^/dy at y=0) and total (oj/^) production cross sections in proton- 
induced interactions. The pA symbol in the second column indicates that the cross section value is 
obtained by fitting different target materials. 

probably in some measurements the systematic uncertainties coming from triggering effects or luminosity 
determination have been underestimated. 

The ifi(2S) cross sections are usually estimated on the same data sample of the crj/^, but are poorer 
in statistics (by a factor around 60). Therefore they are less sensitive to the accurate determination of 
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Experiment 


Reaction 


(GeV) 


(nb/nucleon) 


(J i)(2S)/ (J J/i! 


E331 


28 




pC 


20.6 


15.4 ± 9.1 


0.060 ± 0.035 


E444 


29 




pC 


20.6 


22.8 ± 13.5 


0.137 ± 0.079 


E705 


31 




pLi 


23.8 


42.5 ± 9.0 


0.159 ± 0.029 


E288 


33 




pBe 


27.4 


28.9 ± 11.3 


0.141 ± 0.042 


NA38/51 35, 36 


pA 


29.1 


39.3 ± 9.6 


0.135 ± 0.015 


NA50 23 


pA 


29.1 


47.1 ± 10.9 


0.145 ± 0.017 


E771 


41 




pSi 


38.8 


46.3 ± 5.7 


0.139 ± 0.020 


E789 42 


pAu 


38.8 


66.1 ± 14.1 


0.202 ± 0.055 



Table 3: Updated total production cross sections of ip{2S) an d ratios of if)(2S) to J/tp cross sections in 
proton-induced interactions. 

the systematic effects. Within the quoted uncertainties, there is a good internal consistency among the 
different measurements. 

In the cross section ratios, the systematic effects on luminosity or trigger mainly cancel out and 
the final uncertainty is usually dominated by the ip(2S) statistics. The measurements are all compatible. 

Since the measurements on the polarization are all compatible at 3a with the assumption of no 
polarization ( J401 1411 147j ) . possible effects of the polarization in the J/tjj or ip(2S) decays are always 
neglected and in the evaluation of the cross sections the polarization is assumed to be negligible for both 
charmonium states. It should be noted that if any polarization is present (expecially for ip(2S) where 
a precise information is lacking) some results might be affected by a systematic uncertainty larger than 
that quoted. 

4 Fit results 

The cross sections values obtained from NRQCD calculations at NLO, Eqs. Q and ©, have been used 
to fit the experimental results summarized in Tabs. [21 El 

The ratio of the two color octet matrix elements for J /ip and ip{2S) production to the Tevatron ones 
(\j/ip and Xj/j(2S) respectively) are used to fit the theoretical predictions to the experimental values. 

Since the cross section depends on the product of the matrix elements and the PDF's as shown in 
Eq. (J2), a change in the PDF can influence the result of the fit on the two color octet matrix elements. 
To estimate such effect we used two different PDF sets, the CTEQ6m [IE] and the MRST2002nlo [T9] . 

The data considered in the fit (Tabs. [2] and [3]) include 21 results on the total J/tjj cross section 
(cms-energy range: [6.7 : 200] GeV), three ^(25) cross sections (NA50 [37], E771 [H] and E789 02]; 
cms-energy range: [29.1 : 38.8] GeV) and five cross section ratios R^= &i/ J (2S)/ cr J/ip (cms-energy range: 
[20.6 : 29.1] GeV), for a total of 29 experimental results. This set of data shows a good overall consistency, 
except for the few measurements of cfj/^ already mentioned in the previous section. 

As a first step, the fit has been performed on all the 29 experimental results using the MRST2002nlo 
and the CTEQ6m PDFs. Since the charmonium production from singlet states strongly depends on 
the factorization {hf) and renormalization (hr) scales, we opted for determining the proper scale as an 
additional free parameter of the fit. The resulting optimal scale values are \x = 1.5/Xq for MRST2002nlo 
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and fx = 2.6/io for CTEQ6m, where /i = [if = fiR and fio = 2m c with m c =1.5 GeV. The fit results 
obtained with these scale values are reported in Tab. [I] 





MRST2002nlo 


CTEQ6m 


X 2 /d.o.f. 


1-5/Uo 
114/27 
0.089 ±0.013 
0.061 ±0.012 


2.6/io 
170/27 
0.211 ±0.027 
0.112 ±0.017 



Table 4: Results of the fit performed using the optimal scale factors for the MRST2002nlo and the 
CTEQ6m PDFs. 

The x 2 of the fit is poor for both fits, since, as previously discussed, the measurements of aj/^ are 
not internally compatible. A reasonable explanation is that in some cases the experimental uncertainty 
might have been underestimated. With this assumption, we increased the uncertainties of the fit results 
by a scaling factor (s=2.05 and s=2.51 respectively), obtained following the PDG prescriptions [21] . 
The ratios of the matrix elements (Aj/^,, A^s)) provided by the fits differ by about a factor two mainly 
because of the different scales used. Nevertheless, in both fits the color octet matrix elements are found 
to be much smaller (1/5 — 1/10) than those extracted at the Tevatron ^7j. Given that the CTEQ6m 
PDF shows a worse adaptation to the data, we have decided to choose the MRST2002nlo PDF fit 
as our baseline fit, while the results from the CTEQ6m PDF have been considered in the systematic 
uncertainties determination. 

The fit results as a function of the cms-energy are displayed in Fig. [2]for the MRST2002nlo PDF. In 
the top plot the J ftp cross section is shown. In the bottom-left plot the ip(2S) cross section is shown, 
while the (J^(2S)/ a J/ip ratio is displayed in the bottom-right plot. The open circles in the two bottom 
plots represent the results calculated from the published papers and not used in the fit (see Tab. 01). The 
dot-dashed line indicates the NRQCD predictions without any octet contribution in the charmonium 
production, while the solid line is the result of the fit. The two dashed lines are obtained by independently 
varying the scales fip an d (J>r between hq and 4 (jlq. 

As one can see the curves obtained from the fit can well reproduce the general trend of the J/ ip and 
ip(2S) cross sections as a function of the cms-energy. In order to fit the data the production from singlet 
states seems to be not sufficient and contributions from octet states must be considered. However, the 
fit shows that the amount of octet production is small compared to the results from Tevatron data 
and is in agreement the results on charmonium production obtained at HERA 5 , where the octet 
contributions were negligible. 

The stability of the results with respect to our baseline fit and the systematic uncertainties have been 
checked with detailed studies, mainly by changing the selection of the measurements to be fitted and the 
PDF. First we excluded the experimental results which show a bad partial x 2 (Xp > 9) an d consequently 
are not compatible with the NRQCD calculations. Following these prescriptions the four J/ip cross 
section measurements from the E331 [28] . E705 [SJ, E595 [34j and E771 [41j experiments were left out, 
obtaining a much improved x 2 {x 2 /dof = 37/23), with only a few percent variation on Xjm and A^S)- 
Similar results have been obtained by excluding the four clearly incompatible measurements of crj/^ at 
20.6, 23.8 and 24.3 GeV (E331 [2E], E444 [22], E705 [HI] and UA6 32 J discussed in the previous section. 
More fits have been performed by changing the cms-energy range of the measurements and by using 
different target selections with the MRST2002nlo and the CTEQ6m PDF. Again the results showed a 
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Figure 2: Fit results as a function of the cms-energy for the J/tp cross section (top), the ip(2S) cross 
section (bottom-left) and the (?ip(2S)l (J J/ip ratio (bottom-right). The open circles in the two bottom plots 
represent the results calculated from the published papers which are not used in the fit. 



very good stability on the cross sections; the only significant change observed is, as shown in Tab. 
the dependence of the Xj/^p and \<p(2S) parameters on the actual scales (/ip and /xr) and PDF used. No 
significant variation of the results is observed when we include in the fit, in addition to ojj^ data, only 
a i>(2S) or onr y a i/j(2S)/°'J/ip measurements. In the quoted uncertainties, the correlations between different 
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Figure 3: The J/tp differential cross section da p N /dy at y=0, as a function of the cms-energy. The data 
shown here are not used in the fit. The theoretical prediction with its uncertainty corresponds to the fit 
to the total cross sections, as shown in Fig. El 



experimental results, due to the dilepton decay branching ratios and to the luminosity and efficiencies 
determinations for different measures performed by the same experiment, have been found to have a 
small impact. 

Once the stability of the fit procedure and its precision have been established, one can obtain predic- 
tions for other quantities. As a first application, in Fig. Owe show the comparison between the differential 
cross sections da p ]\r/dy at y=0, as a function of the cms-energy, and the theoretical predictions at NLO. 
The data shown in the plot have not been used in the fit and therefore this result is an important check 
of the overall consistency of theory and experimental data. As a further application of our analysis, we 
present the predictions for the J/tp and tp(2S) cross sections and their ratio at ^/s = 41.6 GeV, where 
HERA-i? [13] has obtained the most recent result on the J/tp cross section from a Minimum Bias data 
sample. The values obtained from the fit are the following: 

ajj^ = (502 ± 44) nb/nucleon , (6) 

a i>{2S) = (65 ±11) nb/nucleon , (7) 

Jfy = (0.130 ±0.019) , (8) 

where the quoted errors include the uncertainties due to the fit, to the data selection and to the PDF. 
The systematic uncertainties due to the PDF have been obtained following the so called "Les Houches 
Accord" prescriptions as discussed in Ref. |48j . 
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5 Conclusions 



In this paper we have collected all available data on J/ip and tp(2S) production in hadron collisions 
(with the exception of data obtained in pion-nucleus collisions) and updated them in the light of the 
most precise determinations of nuclear effects. We have then presented their analysis in the context 
of NRQCD, using NLO predictions for the short-distance cross sections and fitting the color-octet non- 
perturbative matrix elements. In order to ease the comparison with the available determinations, we have 
chosen the values extracted at the Tevatron as our reference. We find sizeable systematic uncertainties 
associated both to the experimental data, which sometimes are marginally consistent among themselves, 
and to the fixed-order nature of theoretical predictions. Nevertheless, our results clearly indicate that 
the amount of color-octet contributions needed to explain fixed-target data is only about 10% of that 
fitted at the Tevatron. One can certainly argue that part of this discrepancy might be associated to 
the fact that the Tevatron analysis is based on LO calculations only. On the other hand, the difference 
is too large to be resolved by the inclusion of higher-order corrections. In addition, it is plausible to 
speculate that once the NLO corrections to the color-singlet production were computed and included 
in the analysis, there would be very little room left for color-octet contributions to fit the fixed-target 
experiments, in close analogy to what happens in photoproduction and in agreement with the results 
in Ref. [IHl where a different approach is used. 

Acknowledgments 

We would like to thank Mike Medinnis and many colleagues of the HERA-£> Collaboration for the 
support during this analysis and for many stimulating discussions. 

References 



[i] 

[2] 
[3] 
[4] 

[5] 
[6] 
[T 



For a review, see N. Brambilla et al., "Heavy quarkonium physics", hep-ph/0412158 . 
E. L. Berger and D. L. Jones, Phys. Rev. D 23 1521 (1981). 
R. Baier and R. Ruckl, Phys. Lett. B 102, 364 (1981). 

G. T. Bodwin, E. Braaten, and G. P. Lepage, Phys. Rev. D 51, 1125 (1995) [Erratum-ibid. D 55, 
5853 (1997)]. 

M. Kramer, Prog. Part. Nucl. Phys. 47, 141 (2001). 
CDF collaboration, Note 05-04-28. 



M. Beneke and I. Z. R othstein, Phys. Rev. D 54, 2005 (1996) [Erratum-ibid. D 54, 7082 (1996)] 
|arXiv:hep-ph/9603400] . 

F. Maltoni, |hep-ph/0007003] . 

G. P. Lepage, L. Magnea, C. Nakhleh, U. Magnea and K. Hornbostel, Phys. Rev. D 46, 4052 (1992). 
[10] G. C. Nayak, J. W. Qiu and G. Sterman, |hep-ph/050902l"| . 



[11] F. Maltoni, M. L. Mangano and A. Petrelli, Nucl. Phys. B 519, 361 (1998). 



11 



[12] A. Petrelli, M. Cacciari, M. Greco, F. Maltoni, and M. L. Mangano, Nucl. Phys. B514, 245 (1998). 

[13] M. Kramer, Nucl. Phys. B 459, 3 (1996). 

[14] W. Buchmiiller and S. H. Tye, Phys. Rev. D24, 132 (1981). 

[15] E. J. Eichten and C. Quigg, Phys. Rev. D52, 1726 (1995). 

[16] F. Abe et al, (CDF Collaboration), Phys. Rev. Lett. 79, 572 (1997). 

[17] P. Nason et al, hep -ph/00031421 . 

[18] H.L. Lai et al, Eur. Phys. J. 12, 375 (2000). 

[19] A.D. Martin et al, Eur. Phys. J. 28, 455 (2003). 

[20] G. T. Bodwin, E. Braaten and J. Lee, Phys. Rev. D 72 (2005) 014004. 

[21] S. Eidelmann et al, Review of Particle Physics, Phys. Lett. B 592, 1 (2004). 

[22] M.J. Leitch et al, (E866 collaboration), Phys. Rev. Lett. 84, 3256 (2000). 

[23] A. Bamberger et al, Nucl. Phys. B 134, 1 (1978). 

[24] M.J. Corden et al (WA39 collaboration), Phys. Lett. 98B, 220 (1981). 

[25] Yu.M. Antipov et al, Phys. Lett. 60B, 309 (1976). 

[26] K.J. Anderson et al, Phys. Rev. Lett. 36, 237 (1976), Phys. Rev. Lett. 37, 799 (1976). 

[27] J. Badier et al, (NA3 collaboration), Z. Phys. C 20, 101 (1983). 

[28] J.G. Branson et al, (E331 collaboration), Phys. Rev. Lett. 38, 1331 (1977). 

[29] K.J. Anderson et al, (E444 collaboration), Phys. Rev. Lett. 42, 944 (1979). 

[30] J.H. Cobb et al, Phys. Lett. 68B, 101 (1977). 

[31] L. Antoniazzi et al, (E705 collaboration), Phys. Rev. D 46, 4828 (1992). 
[32] C. Morel et al, (UA6 collaboration), Phys. Lett. B 252, 505 (1990). 
[33] H.D. Snyder et al, (E288 collaboration), Phys. Rev. Lett. 36, 1415 (1976). 
[34] E.J. Siskind et al, (E595 collaboration), Phys. Rev. D 21, 628 (1980). 
[35] M.C. Abreu et al, (NA51 collaboration), Phys. Lett. B 438, 35 (1998). 
[36] M.C. Abreu et al, (NA38 collaboration), Phys. Lett. B 444, 516 (1998). 
[37] B. Alessandro et al, (NA50 collaboration), Eur. Phys. J. C 33, 31 (2004). 
[38] C. Kourkoumelis et al, Phys. Lett. 91B, 481 (1980). 
[39] A.G. Clark et al, Nucl. Phys. B 142, 29 (1978). 



12 



[40] A. Gribushin et al, (E672/706 collaboration), Phys. Rev. D 62, 012001 (2000). 

[41] T. Alexopoulos et al, (E771 collaboration), Phys. Rev. D 55, 3927 (1997). 

[42] M.H. Schub et al, (E789 collaboration), Phys. Rev. D 52, 1307 (1995). 

[43] I. Abt et al, (HERA-B collaboration), submitted to Phys. Lett. B., hep-ex/0 5 12029 . 

[44] E. Nagy et al, Phys. Lett. 60B, 96 (1975). 

[45] E. Amaldi et al, Lett. Nuo. Cim. 19, 152 (1977). 

[46] S.S. Adler et al, (PHENIX collaboration), Phys. Rev. Lett. 92, 051802 (2004). 

[47] T.H. Chang et a/.,(E866 collaboration), Phys. Rev. Lett. 91, 211801 (2001); T.H. Chang, PhD 
Thesis, |hep-ex/0012"034l 

[48] J. Pumplin et al, JHEP 0207, 012 (2002). 

[49] V. A. Khoze et al, Eur. Phys. J. C39, 163, (2005). 



13 



